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We review some properties of the stationary states of the Fokker - Planck equation for TV 
interacting particles within a mean held approximation, which yields a non-linear integrodifferential 
equation for the particle density. Analytical results show that for attractive long range potentials 
the steady state is always a precipitate containing one cluster of small size. For arbitrary potential, 
linear stability analysis allows to state the conditions under which the uniform equilibrium state is 
' unstable against small perturbations and, via the Einstein relation, to define a critical temperature 

T c separating two phases, uniform and precipitate. The corresponding phase diagram turns out 
to be strongly dependent on the pair-potential. In addition, numerical calculations reveal that the 
transition is hysteretic. We finally discuss the dynamics of relaxation for the uniform state suddenly 
cooled below T c . 
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I. INTRODUCTION 
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The dynamics of brownian particles is a subject of great interest in statistical physics, especially when the 
particles interact through potentials that can be short or long ranged, from the simple hard-core to coulombian 
interactions. The physics of surfaces, more precisely the motion of adatoms on substrates, provides an experimental 
I | realization of this problem. Interactions can lead to collective phenomena, and, generally speaking, to patterns in 
space. Pioneering works have re,vfaled experimental and computational evidences of a phase transition in the case 
of oxygen adsorbed on tungstentra. This phase transition and the dynamics of such systems (e. g. the modification 
of the diffusion constant I? of a tracer, time auto-correlation fonction of the on-site density, etc. . . ) where analyzed 
i__T by Tringides et alB and reviewed by GomerQ- Numerical studies for a two-dimensional lattice gas have also been 
performed in the case of a contact interactiorJjQ. On a different scale^the hydrodynamics of interacting brownian 
particles have been studied for short-ranged or screened interactions inErtl In a slighly different context, especially 
when the-eonsidered system is an open one, many papers have been devoted to statistical models describing chemical 
reactions!!!] (and references therein, especially the review by Zhdano and KasemcO). Ear details on 2-dimensional 
lattice gas models (and their subsequent approximations), the review by Kehr and Binderlij is highly recommendable. 

In the present paper, we use the framework of the many-body Fokker-Planck equation (FPE), written for N 
interacting particles - as such, this is a continuous space model. The ordinary FPE describes the diffusive motion of 
, brownian .particles under the assumption of slow diffusion. It can be derived by using the Kramers-Moyal systematic 
expansion^ and arises when one assumes that all the moments of the increments of the stochastic variable of order 
> 3 are proportional to Ai r , r > 1, where At is the time increment. As such, FPE is a conservation equation for a 
probability density P, which can always be written in the form: 



a t p = -divj, (i) 

where J is the probability current. As an example, for a single particle with position x moving in the static external 
potential ^(x), the current is the sum of the drift term — /j,Pd x V and of the diffusion current —Dd x P; in such a case, 
' the equation writes: 

■ 5 * P ( X ' *) = d * i D 9 * P ( x ' f ) + A* P ( x ' *)3xV(x)] , (2) 

where /i is the mobility and D the diffusion constant. Eq. (g) gives the probability density for the position of a 
particle obeying the Langevin equation in its viscuous limit, an initial distribution being given. 

With several particles, the interesting case occurs when they interact through a given internal force field deriving 
from a potential V , opening the possibility of competing effects. When V is purely repulsive, no interesting effect is 
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expected: in infinite space, one can guess that the equilibrium state is the uniform one, all probability densities being 
constant in space. In the opposite case, when the particles attract each other, a competition between diffusion and 
interaction takes place, which can, at least in principle, produce patterns or structures in space. Obviously enough, 
the possibility of the latter depends on the features of the interaction potential, namely its strength and its range, 
and possibly of the dimensionality. The competition between drift and diffusion is measured by the ratio D / fj,; as a 
consequence, patterns can be expected at low temperatures, i. e. when the drift term dominates diffusion. If there 
exists a definite value D c of D which separates two distinct stationary solutions, the Einstein relation D/fj, = IcbT 
allows to identify a critical temperature T c . 

The purpose of this paper is to put forward a few results concerning the equilibrium state of FPE for interacting 
particles, obtained within a mean field approximation. The paper is organized as follows. After setting the basic 
equations relevant to our purpose, we first focus on some potentials allowing an exact treatment of the mean-field 
equations. In a second step, we discuss the linear stability of the uniform equilibrium state, which is always a solution 
of the problem. It is seen that unstabilities can indeed occur, and the conditions for that are given, yielding the 
expression of the critical temperature T c . Eventually, we give the far-from-equilibrium dynamics of the uniform state 
suddenly cooled below T c . 

For N identical interacting particles, the potential is noted V(xi, x 2 , . . . xjv) and the generalization of (||) 

writes: 



<9 f P(xi, . . . , xjv, t) = ^9 Xi [Dd Xi P{xx, . . .,x N ,t) +fiP(xi, x N , t) d Xi V(xi, 



(3) 



in the following the potential is assumed to be the sum of N(N — 1) even two-body terms: 

V(xi, x 2 , . . . Xjv) = - ^2 v ( x i ~ x j) > v ( x ) = v (~ x ) ■ 



(4) 



Obviously, the solution of (|J) - if known -, is of little physical interest, since one is usually interested in the 
one-particle density P^ and the pair-density function P^ (reduced densities of order 1 and 2), defined as: 

P (1) (x, t) = J dx 2 . . .d XA rP(x, x 2 , . . . , xat, t) , P( 2 >(x, x', t) = J dx 3 . . .dxArP(x, x', x 3 , . . . , Xat, t) . (5) 

Due to the many-body interactions, reduced densities obey a hierarchy of the BBGKY type (see e. g. McQuarrieEl) ; 
in the present context, the first equation of this hierarchy writes: 



D <9 X P (1) (x, t) + (N - 1) n J dx' PW (x, x', t) d x v{x - x') 



(6) 



Solving this hierarchy is usually impossible; the simplest approximation is of the mean-field type, in which one imposes 
the form: 



P(xi, x 2 , ... xat, <) = JJ p(x 4 , t) . 

i=l 

From @ , the mean-field one-particle probability density p must obey the following equation: 



d t p(x, t) = <9 X 



D 9 x p(x, t) + (N-l)n p(x, t) J dx' p(x', t) 9 xW (x - x') 



(7) 



(8) 



Thus, for N 3> 1, the particle density n(x, t) — Np(x, t) is the solution of the integrodifferential non-linear equation: 



9 t n(x, t) = <9 X 



D d x n(x, t) + fi n(x, t) d x J dx' n(x', t) w(x — x') 



(9) 



Clearly, even in this mean-field approximation, finding the solution is by far not a trivial question. 
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II. SOME EXACT STATIONARY SOLUTIONS OF THE ONE-DIMENSIONAL MEAN-FIELD FPE 



In one dimension, all the stationary solutions of (||) give a vanishing current J and are the solutions of: 



d d /*+°° 

D — n(x) + n n(x) — / dx'n(x')v(x 
dx dx J^^ 



= . (10) 



Because n also appears in the integral, the stationary mean-field solution is not connected in an obvious way to the 
Boltzmann distribution built with v(x). One could naively believe that, since in the mean-field treatment each particle 
interacts with N — 1 ~ N other particles through the potential v(x), its equilibrium distribution is oc e ~ N v ( x ) _ This 
turns out to be wrong in general, except for the harmonic potential. Indeed, the one-particle density is obtained by 
integrating over the coordinates of all the other N —1 particles and there is no reason ensuring that the bare two-body 
interaction v(x) should spontaneously appear in the Boltzmann way in the one-particle density, even in a mean- field 
approach. It will be seen that the ordinary Boltzmann factor is recovered only in certain limits (see below). Also note 
that, since there is no external force field, the equilibrium states - assumed to be independent of the initial condition 
which naturally implies privileged points - are defined up to an arbitrary translation in space. Otherwise stated, if 
n(x) is a solution of ( |To| ) defined for all x between ±oo, then n Xo (x) = n(x — xo), xo arbitrary, is also a solution. 
This degeneracy is discarded by the use of symmetric boundary conditions (see below) or by having in mind that the 
displayed equilibrium states arise from the initial condition p{x, t = 0) = 8{x). 

For some potentials, the solution of the equation ( |l0| ) can be found in closed form. In the following section, we 
give a few examples and briefly analyze the corresponding solutions. 



A. The Coulomb Potential 

We first choose: 

v(x) = v M , (v Q > 0) , (11) 

which, for d = 1, mimics the Coulomb potential in the sense that v(x) satisfies the Poisson equation. This is clearly a 
long-ranged attractive potential: the force exerted on a given particle is constant in space and is equal to — (i>o/£) sgnx. 

Let us first assume that the particles are confined in the interval [-L/2, +L/2]. Introducing the integrated 
density q(x): 

q(x) = [ n(x') dx' , (12) 



it is readily seen from ( |To| ) and (|Tl]) that q(x) obeys the following differential equation ((3 1 = fee? 1 ): 

q"(x) + ^q'(x)[2q(x)-N] = . (13) 

Using the boundary conditions q(—L/2) — and q(+L/2) — N, a little algebra yields the properly normalized density 
in the limit L — » +oo: 



2 A cosh 2 (x/ A) ' N(3v 



= TTi TTT . A = . (14) 



This is a peaked distribution, which tends towards N8(x) in the limit N — > +oo. Due to the infinite range of the 
potential, the stationary state is for any D and large N a cluster of very small shape as compared to £ when (5vq ~ 1. 
Note that p(x) is flat at x = but decreases approximately like an exponential: 

n(x) ~ r X N 2 /3v e- N ^ x ^ = r 1 N 2 (3v e ~ N ^ if x » -|- • (15) 

Npva 

This shows that the Boltzmann behaviour is recovered in the wings of the distribution only. 

Note that it is not necessary to introduce a finite interval of length L and subsequently to take the limit L — > oo. 
Yet, this procedure allows to discuss the invariance by translation mentionned above. Indeed, taking the boundaries 
at ±L/2 forces the solution n(x) to be even, a symmetry which is conserved when the limit L — > oo is performed 
afterwards. On the other hand, by assuming infinite space at the beginning, the calculation yields the same solution 
( |l4| ) as well as all the translated functions n Xo {x) — n(x — xq) with xq arbitrary. It can be checked that all the 
functions n Xa (x) indeed satisfy (|To|). 
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B. The Harmonic Potential 



We now take v(x) = v x 2 / (2t; 2 ) , vq > 0, another example of long-ranged attractive potential. In this case, ( |l(i| ) 

writes: 



Pv 

For even n(x) = n(—x), this simplifies to: 



n'(x)+ n(x) / {x-x') n(x')dx' = . (16) 



— OO 



n'(x) + N^-xn(x) = , (17) 



which gives the solution for the density: 



n{x) = _^ e -* 2 /(2A 2 ) A = ^jL= . (18) 



Note that the expression (|l8] ) is simply of the form oc e~ N l 3v ( x )^ which is the Boltzmann distribution for a single 
particle elastically bound with N others. The fact that this is true for all TV and a; is clearly characteristic of the 
harmonic potential. Again, it is readily checked that all the functions n(x — xo) satisfy (E6h with xq arbitrary. 



C. Polynomial Potentials 



The harmonic potential treated above immediately gives the clue for solving the same problem with any 
potential of the polynomial type: 



v(x) — c r x r . 

r>l 

With this kind of potential, ([l0]) assumes the form: 

/+oo 
{x-x'Y" 1 n(x')dx' = . 

By rearranging terms in the integral, this can be rewritten as: 

n'(x) + (3 n(x) ^ lr(x) f r = , 

where the "f r (x) are definite polynomials and where the quantities f r are the moments of n(x): 



fr = 



x r n(x) dx 



(19) 



(20) 



(21) 



(22) 



Now, ( |2l| ) can be formally integrated, to give a function p(x) containing the parameters f r : 

n{x) = e -PE r> _ ^^ f (23 ) 

the T r being definite functions depending on v(x). By reporting this expression in (^2|), one can write as many 
equations as necessary to find the f r and eventually obtain the explicit expression of n(x). Clearly, the latter is not 
a priori of the form oc q- n v ( x ) . 

As an example, let us consider the quartic potential: 



v(x) = Vq 



X 

9 I T 



(24) 
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This potential is purely attractive if g > 0; otherwise, it is repulsive for x between g/2] and attractive elsewhere. 

Inserting this potential in ( ^o|) and integrating, one finds: 

n{x) = Ce-" , °« 9 ^° W2 > l! /« 4+/ » l4 / $4 l . (25) 
/o is equal to N, whereas the unknown quantities C and f% can be derived from the two equations: 

p + co p + co 

2C / e -Pv [{gN+f 2 /e)X 2 +NX i ] AX = N ^ 2CS 2 / ^2 & -0v o [(gN + f 2 /e)X 2 +N X l \ dJ{ = ^ ^gj 



Note that the expression (pia) is not c>c e N l 3v ( x ) ; except if /a is negligible as compared to N^ 2 (see below). The above 
integrals can be expressedEJ with the Weber functions T> v : 

^C(2N{3v )- 1/4 e A2 /*V_ 1/2 (A) = N , ^ C^ 2 (2iV/?« )- 3 / 4 e^/ 4 2?_ 3/2 (A) = / 2 , (27) 
where the constant A is: 

A= ^( 5 ^ + / 2 /£ 2 ) 2 ■ (28) 

For arbitrary N, these equations can be numerically solved to provide the two parameters C and / 2 . On the other 
hand, when N is very large, one can find directly asymptotic expressions of the integrals appearing in (^). By doing 
so and coming back to the mean-field density n(x), one gets according to the sign of g: 

n{x) ^ NC 1 (j^f-) V2 e- NMx) (9 > 0) , (29) 



1 /2 

n(x) ~ NC 1 (^=^) e- w ^ 2 / 4 e-^) (g < 0) . (30) 

This shows that the Boltzmann distribution is recovered only in the limit N 3> 1 . For positive g, one has a single 
cluster with a width ~ (NPvq)^ 1 ^ 2 £. For negative g, two peaks arise, both having the latter width, and separated 
by x 

Note that, in any case with (3v$ ~ 1, the mean-field stationary state for large AT is a compact cluster with a 
very small size as compared to the characteristic length £ of the potential. 



III. LINEAR STABILITY ANALYSIS OF THE UNIFORM STATE 



In the previous section, we displayed several potentials allowing an exact explicit expression for the non-uniform 
equilibrium state with a vanishing current. On the other hand, ( |To| ) always trivially has the uniform state as solution. 
All this means that several solutions can exist and the question arises to settle which of them can indeed be realized 
for a given ratio D/fi. The aim of this section is to analyze the linear stability of the uniform solution no of the 
mean-field equation ( |Io| ) in any dimension. Setting n(x, t) = no + Sn(x, t) and discarding all terms of order greater 
than one, one obtains: 

d t Sn(x, t) = Ddl Sn(x, t) + fin d x J dx' Sn(x', t) <9 X u(x - x') . (31) 
This linear equation can now be analyzed by introducing the Fourier transforms: 

p(q,t) = / dxe- iqx <5n(x, t) , u(q) = / dxe~ iqx i-(x) , (32) 



and it is readily seen that all the eigenmodes of (|T]) are of the form e "( q '*, with w(q) given by the dispersion relation: 

w(q) = [D + u-n i,{q)} q 2 = £> off (q)q 2 . (33) 
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D c g plays the role of an effective q-dependent diffusion constant; it can be said that, due to interactions, Fick's law 
is no more a local law. Setting: 

v(x) = Vo4>(x/£) , 0(0) > , x = |x| , q = |q| , (34) 

the dispersion law (|33| ) in d dimensions can be rewritten as: 

uj(q) = D [l+pv n ^UM)] ? , (35) 
where Ud is a dimensionless function. More precisely, one has (Jo is the ordinary Bessel function): 

r + co r+oo 

C/i(k) = 2 / (f>(X) cos(nX) dX , U 2 (n) = 2tt / X 4>{X) J {nX) dX , (36) 



U 3 (k) = — / X(f>(X) sm(nX)dX . (37) 
« 'o 

Eq. ( |35| ) shows that the uniform solution is unstable against a deformation with a wavevector q satisfying: 

k B T + n C d v UM) < ■ (38) 

Clearly, such a condition cannot be satisfied for a purely repulsive potential, since then both vq and Ud are positive 
quantities^]: for repulsive forces, the uniform state is always stable, as anticipated on physical grounds in section |[ 

When v is attractive at short distances, vo is negative. The possibility of instability then depends only on the 
precise behaviour of the dimensionless Fourier transform Ud- For very short ranged potentials, Ud{n) is expected to 
be positive for all n. The instability condition can then be rewritten in a more transparent form : 

k B T < n ^ d \v Q \U d (qO , v Q <0 . (39) 

If, in addition, Ud is a monotonous decreasing function, instability occurs at large wavelengths provided the temper- 
ature is low enough. More precisely, provided that Ud(0) is finite (short-range potential), one can define a critical 
temperature T c : 

k B T c = n o C d \v o \U d (0) . (40) 

For T < T c , there exists a finite interval [0, q m &x] in which the uniform state in unstable. g max is a function of the 
density no and the temperature. Assuming that Ud(0) is of order unity, ( |4C| ) means that the critical temperature is 
such that the thermal energy is of the order of the effective interaction energy. 

The existence of an upper g max , below which the uniform solution is unstable is easily understood on physical 
grounds. For a small disturbance having a large wavelength compared to the range of the potential, particles in excess 
tend to attract each other more firmly, increasing holes in the density. On the other hand, when the wavelength is 
quite small, the inhomogenities arc easily washed out since holes and particles all interact. 

As an example, let us consider the d-dimcnsional gaussian attractive potential: 



In this case, one has: 



and q max is given by: 



v (x) = v e- x /(2 « > , vo < . (41) 
k B T c = (2ir) d / 2 n t d \v \ (42) 



1/2 



= rM 21n ^) . T<T C (43) 



1 When 4>(X) is a positive function, Ud is bounded below by a positive (possibly diverging) integral. 
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which entails that: 

gmax ~ (T C -T)V 2 , T <T C . (44) 

These last results are depicted in the left part of fig. [l| The same conclusions still qualitatively hold if a repulsive 
core is added (e. g. 6-12 Lennard-Jones or Morse potential). This only affects the large-g behaviour of v and does 
not change the results. 




FIG. 1. Schematic representation of the instability domains for the gaussian potential (^) (left, below <jr ma x) and for the 
square potential (|^|) (right, pockets). 

When the potential has a sharp cut-off, its Fourier transform has an oscillatory behaviour at relatively small 
wave numbers and an interesting phenomenon occurs. To be specific, let us take the following square potential: 



v(x) 



vo if M < f/2 
otherwise 



v < 



(45) 



We here find k^T c — noS,\vo\- Now, let us note the abscissa? of the maxima of (smX)/X and Sk = (sinXfc)/Xfc 
(Xo — 0, so — !)• Then, for s\T c < T < T Cl a single instability interval arises, [(?min, l = 0, <7max, l 1. When T 
decreases, another interval [g m m,2, <7max, 2] is found when S2T C < T < s\T c , and so on. The situation is depicted 
on fig. |^. Thus, for a short-range potential with a sharp cut-off, several disjoint instability intervals are successively 
obtained when the temperature is decreased. 

If Ud{0) is infinite, T c is formally also infinite; this means that the uniform solution is unstable at small q 
for any temperature (see fig. |[ left). For instance, for the one-dimensional Coulomb potential (11), one has (after 
regularization) J7i(«) = 2/k 2 , so that: 



r 1 



k B T 



(46) 



The same holds true for higher dimensions, since Ud{n) is for any d proportional to k 





FIG. 2. Schematic representation of the instability domain for the Coulomb potential ( |ll| ) (left, below q„ 
Morse potential (^?j) (right, between g m i n and g ma x). 



and the inverted 



In the preceding examples, the first instability interval, when it exists, arises around q = 0. This is, as explained 
above, a physical consequence of the fact that the potential is purely attractive (the fact that T c is finite or not is a 
consequence of the range only). Another interesting situation is when v is attractive at short distance and repulsive 
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at long distance; this means that two particles are bound by a potential barrier but can dissociate when one of them 
is given enough energy. As an example, let us consider the inverted 3d Morse potential: 



v(x) 



-2(r-r )/€ _ 20-^-^^ 



v < 



(47) 



In such a case, T c is again finite, and the instability interval for q is of the form [q m in, 9max]; it grows around a finite 
value (which gives U^^q^) its secondary minimum) and is a consequence, at intermediate wavelengths, of the interplay 
between attraction at short distances and repulsion at large ones. This instability domain is schematized in figure |^. 
Quite naturally, the uniform state is stable against long wavelength deformations (due to the long distance repulsive 
behaviour of the potential) and unstable otherwise. 



IV. NUMERICAL STUDY OF THE GAUSSIAN ATTRACTIVE CASE 



From the arguments and results given in the two previous sections, we expect that for a Gaussian attractive 
(an example of short-range potential), two equilibrium solutions exist; the uniform one, which is unstable against long 
wavelength deformations when T < T c , where T c is given by (^2j), and a non- uniform distribution, which we were 
unable to find analytically in a closed form. The aim of this section is to give numerical results for the one-dimensional 
gaussian potential. It is shown below that the predicted transition between these two states is indeed quite sharp 
when T crosses T c and that, in addition, an hysteretic behaviour occurs. 

We choose v(x) = v e~ x ) (vq < 0) which yields v(q) = y/2ir £vq q I 2 . The unstability condition ( |38| ) 
gives the critical temperature (see 



2tt ra £NI/fc 



(48) 



In order to analyze the equilibrium state(s) n(x) arising in the mean-field treatment with the gaussian potential, 
we used two numerical procedures. The first one is a self-consistent iterative procedure for the equilibrium equation 
(|l0|), hereafter called static algorithm. The second one (called dynamical algorithm) numerically solves the one- 
dimensional version of the full equation (^), starting from a given initial condition n(x, t = 0), and gives the non- 
equilibrium evolution of the system. 

With the static algorithm, we start from a trial density function and then iterate according to: 



n {k+1) {x) = Z~ 1 e~ /3 -^ dj/ n{k) ( x Mv-x) 

where Z is the normalisation constant. Convergence is reached when the following quantity: 

f f n ( - k+1 \x)~n < - k \x) \ 2 
J L dx { n(k+D(x) ) 

becomes much less than unity. Strict boundary conditions have been used, namely j(0,t) = j(l,t) = 



(49) 



(50) 



n(x) 



n(x) 




FIG. 3. Equilibrium states for the gaussian potential (left: T/T c = 0.9, right: T/T c = 1.1); the length of the box is L = 100£ 
and N = 100 (no = 1). On the left curve, the crosses are the results from the static algorithm whereas the solid line is the 
gaussian approximation denned below. £ is the unit length; note the different horizontal and vertical scales. Above the localized 
curve is shown the lenght £. 
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These calculations allow first to check the existence of the sharp transition occurring at T c when the temperature 
is decreased starting with the uniform state; second, they reveal a hysteretic behaviour: when the temperature is 
increased starting from the localized state (one cluster), the transition from the localized state to the uniform one 
occurs at another temperature T* > T c . 

Fig. [3] displays two typical density profiles (uniform vs localized cluster), just above and just below T c obtained 
with this static algorithm. Since the width of the localized state turns out to be a little bit smaller than the range £ 
of the potential, the central part of the cluster can be safely approximated by a gaussian; indeed, when one expands 
the potential near x/£ ~ 0, the harmonic potential is recovered (if instead of a gaussian potential, one chooses 
v(x) = woe - ' 1 '^: the central shape of the localized cluster is 1/ cosh 2 x, as it must (see ([l4j)) since near the center 
this potential is essentially the Coulomb one). Obviously, such an approximation cannot properly represent the tails 
of the equilibrium density. 

In order to characterize the transition, we introduce two quantities, e and A defined as follows: 



1 



N 



e(T/T c ) = — / dx dx' n(x)n(x')v(x — x') 



(51) 



and: 



A 2 (T/T C ) = — J dx x z n(x) 



N 



dx x n{x) 



(52) 



e represents the mean-field interaction energy per particle; A measures the spatial width of the cluster. For the 
uniform state one has e U nif = {N/L) J L dxv(x) ~ N(£/L) vq and A un jf ~ L ^> £; on the contrary, for the localized 
state, eioc ~ Nvo and Ai oc < £. For an infinite system, the ratio eioc/eunifonn vanishes, so that s can play the role 
of an order parameter. 

The hysteresis is obtained as follows. For each given temperature, the starting point is the solution found in 
the previous run at a slightly different temperature. After a first sequence where T/T c is decreased step by step down 
to 0.2, the procedure is reversed, T/T c is increased, the iterative procedure starting again from the localized solution 
obtained for the previous value of T/T c . It is found that the inverse transition does not occur at T/T c = 1, but at a 
higher temperature T* . 

Clearly, for T < T c , several solutions can exist. Since the width of the single peak solution is a bit smaller than 
the range £ of the potential, one can guess that there also exists a solution with two peaks separated by a distance 
much larger that £. This is confirmed by the numerical calculations: starting with a trial function having two separate 
peaks, the calculation converges towards a stable solution having the same features. One thus can claim that, on the 
low-temperature side, localized solutions exist displaying j = 1, 2, 3, . . . peaks. Each of these solutions has its own 
temperature T* ■ . The corresponding typical hysteretic cycles are shown on figure (^) for the "order parameter" e for 
the one-,two-, and three-peak solutions. 

The possibility of the n-peak solutions is confirmed by the dynamic algorithm (it was checked that the latter 
eventually yields the same equilibrium solutions as the static one). This allows to conclude that a given j-peak 
solution is linearly stable for T c < T < T* ^ . Note that the linear stability of the n-peak solutions above T c is not in 
contradiction with the results obtained in section III; all this simply means that, in the region T c < T < T* ■, both 
the uniform state and the localized state(s) are stable against small perturbations. 



100 T _ 
80-- 




FIG. 4. \e — e un if| versus T: hysteresis cycle for the uniform-localized transition. Solid line: 1-peak branch, dashed line: 
2-peak branch, dotted line: 3-peak branch. 







Defining the entropy as S = — J L dxn(x) ]xi[n(x)/no], it is seen that the one-peak solution is the one having 
the lowest free energy F = e — TS. It must be noted that each kind of states has its own relevant physical parameters. 
For the uniform state, the relevant parameters are vqj^uo = N/L and the critical temperature found by the linear 
stability analysis is a function of these parameters. For the localized state, the relevant parameters are vo,£ an d ^V", 
the number of particles. Indeed, for this state, edge-effects become irrelevant when («L: in this sense L becomes 
arbitrary so that no and N are in fact independent variables. 

The physical reason of the hysteresis can be understood as follows. In the uniform state, a given particle 
interacts with a small (microscopic) number of particles (~ no£); on the contrary, in the localized state, a single 
particle interacts with a large number (~ N), macroscopic in the sense that all the particles effectively interact with 
anyone of them (note that this also originates from the fact that the particles are of zero radius). The condition for the 
uniform-to- localized transition is expressed by k^T c ~ ^oCl w ol ( see (p5|))- For the inverse transition, one can expect 
that a somewhat similar condition also holds true, so that T* is an increasing function of N (and is independent of 
no)- It is likely that with pointlikc particles T* goes to infinity with N; on the contrary, with finite-size particles, 
one can figure out that T* saturates at very large N. As a whole, the hysteretic behaviour is a consequence of this 
parameter cross-over between the two kinds of steady states. 

These facts are illustrated by looking at the width A of the cluster in the localized state, displayed in fig. || 
where A is plotted as a function of the relevant variables for the uniform state (as they naturally arise in the linear 



stability analysis of section III). Indeed, choosing N = 50 or N = 100 with fixed no gives two different curves. 

A/5 



FIG. 5. A versus T/T c , for N = 50 (solid line), and N = 100 (dashed line). The two curves display two well separated 
regimes; at low temperature, A oc T 1 ' 2 , whereas at high temperature, A still follows a power law, but with an exponent ~ 5.4. 
Both curves end at their own localized-uniform transition temperature. 

In the localized state, the width A of the cluster displays a cross-over between two regime s. O n the low-T side, 



one finds approximately A cx T 1 / 2 . This behaviour is the same as for the harmonic case (see [I B), which suggests 
that a gaussian approximation is valid in the low-T region. Indeed, one can compare A with the analytic expression 
( [ilf ) obtained for the harmonic case (see figure |6|) . The agreement turns out to be quite good and this allows to claim 
that a gaussian approximation provides a very accurate description of the localized state arising in the low-T region. 
For the second regime (T c < T < T*), the width of the cluster becomes of the same order or even greater than the 
potential range, so that the effective potential cannot be well represented by a harmonic approximation. Indeed, the 
increase of the cluster width with T no more displays a gaussian behaviour, but still follows a power law. 

A/5 



FIG. 6. Comparison between the gaussian approximation (dashed line) and the numerical results (solid line) for the depen- 
dence of the width vs the temperature, in the low-T region (TV = 100). 
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V. DYNAMICS OF A SYSTEM OF PARTICLES INTERACTING THROUGH THE GAUSSIAN 

POTENTIAL 



In the previous mean-field analysis, it was shown that brownian particles interacting through a gaussian po- 
tential can have several steady states when T c < T < T*. The system will converge towards one of those states, 
depending of the initial condition. It can be expected that, under a strong perturbation, the former is able to shift 
from one steady state to another, without changing the two phase space parameters (T/T c and N). 

We here give a few results obtained with the dynamical algorithm for a set of N = 100 particles starting in the 
uniform state at T/T c = 0.5. The length of the box is L = 100£. This represents a case where the system, being in 
its stable uniform state at a given T > T c , is cooled infinitely quickly below its critical temperature. 

We find that the dynamics is a succession of steps, each of them being a metastable state. Each plateau has a 
free energy lower than that of the preceding. The mean energy per particle is plotted in fig. [?| as a function of time. 
In order to be complete, some corresponding density profiles at different times are shown in fig. H. 



£ 
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FIG. 7. Time evolution of the mean energy per particle e, starting with the uniform state sudenly cooled below the critical 
temperature (T/T c = 0.5 and N = 100). Note the horizontal logarithmic scale; the time unit is £/(/mo|uo|)- 

These plots show that the dynamics is a succession of metastable states with a given number of peaks of high 
density, (13-peak state, 7-peak state, 5-peak state, 3-peak state). The lifetime of these states rapidly increases when 
the number of peaks gets lower (note the horizontal logarithmic scale in fig. ^). 

Because of computational time, we were unable to obtain the 1-peak state. From static computations, it is 
believed that this state is expected to be the final one, although there is no evidence that the two algorithms have 
the same attraction basins. Moreover, the 1-peak state is not the only stable state: we have checked that also exist 
stationary 2-peak, and 3-peak states, a fact which seems nearly obvious on physical grounds since when the distance 
between two neighbouring peaks is much larger that the range £ no dynamics can occur. 

From the time-dependent density profiles (fig. |J), we see that "nucleation" always occurs at the boundaries 
of the box. This is easily understood, since a particle near this boundary is pulled from one side only, as opposed 
to a particle located in the bulk. This effect would disappear with periodic boundary conditions, indicating that the 
nucleation can only occur around an inhomogeneity. 

It is also interesting to follow the annihilation of the peaks. As an example, one goes from 13 to 7-peaks via 
the diffusion of the even peaks into the odd ones. Generally speaking, it turns out that the transition towards a 
state having fewer peaks proceeds through the absorption of "weak" peaks into strong ones. The actual value of the 
density between the peaks, even when it is extremely small, turns out to be of first importance in this mechanism of 
coalescence. 
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FIG. 8. Sequence of the density profile at different times, starting from the uniform state cooled below T c . 
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VI. CONCLUSION 



In this paper, we have studied properties of the iV-body Fokker-Planck equation within a mean field aproxima- 
tion. Our aim was to analyse the properties of the steady state - which is a uniform density gas when interactions are 
set to zero - as a function of the features of the interaction potential. The linear stability analysis of the uniform state 
allowed to state the condition for a phase transition, as a result of the competing effects due to attractive interaction 
and to diffusive spreading. In particular, it was shown that the uniform phase is always unstable for non-summablc 
attractive potentials, the critical temperature being infinite. In this latter case, analytical results have been given for 
specific potentials (ld-coulombian, harmonic, polynomials), all yielding an aggregated steady state. 

For an attractive short ranged interaction, the above stability condition was interpreted as a phase transition 
occurring at a finite T c . This transition was more thoroughly studied by numerical computations in a definite 
case (gaussian potential) and turns out to be of the first order and hysteretic, the low temperature state being an 
inhomogeneous state having one or several peaks of high density. 

The far-from-equilibrium dynamics was numerically computed for the gaussian potential, starting from the 
uniform state suddenly cooled below T c . The relaxation proceeds through a sequence of plateaux displaying a given 
decreasing number of peaks, each of the former having a longer and longer lifetime as time goes on. Each of these 
steps can be viewed as a metastable state. 

Obviously enough, the considered problem requires further investigations. One of them is the relevance of mean 
field approach. If one expects that mean field results should be basically valid for long ranged potentials and high space 
dimensionality, their correctness for short range potential has to be checked by using more refined approximations 
from iV-body general methods. 
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